 SUBROUTINE TRANSPOSE_YZ_XZ(F,nxpsi,nxpse,nypsi,nypse,nzpsi,nzpse, &
     &                     Q,imin,imax,jmin,jmax,kmin,kmax,nx,ny,nz, &
     &                     comm3d,nproc,ndim,dims,coords, &
     &                     nxpi,nxpe,nypi,nype,nzpi,nzpe,isize)
 
!This subroutine transposes the YZ array decomposition into
!XZ decomposition

! use MPI
 use MPI_VARIABLES, only : &
 & transpose_yz_xz_indic
 use TIMING_MOD, only : &
 & proc_time

 implicit none

 include 'mpif.h'


!Input variables 
 real(8), intent(in) :: F(nxpsi:nxpse,nypsi:nypse,nzpsi:nzpse)
 integer, intent(in) :: nxpi(0:nproc-1) ! The boundaries of the portion of array being exchanged
 integer, intent(in) :: nxpe(0:nproc-1)
 integer, intent(in) :: nypi(0:nproc-1)
 integer, intent(in) :: nype(0:nproc-1)
 integer, intent(in) :: nzpi(0:nproc-1)
 integer, intent(in) :: nzpe(0:nproc-1)
 integer, intent(in) :: nxpsi,nxpse,nypsi,nypse,nzpsi,nzpse
 integer, intent(in) :: imin, imax
 integer, intent(in) :: jmin, jmax
 integer, intent(in) :: kmin, kmax
 integer, intent(in) :: nx, ny, nz
 integer, intent(in) :: ndim, nproc
 integer, intent(in) :: comm3d
 integer, intent(in) :: isize(0:nproc-1) ! Corresponds to basic x-decomposition of the exchanged array     

!Output variables
 real(8), intent(out) :: Q(imin:imax,jmin:jmax,kmin:kmax)

!Local variables 
 integer :: i, j, k
 integer :: ki, icoord, max_maxbuf
 integer, intent(in) :: coords(1:3)
 integer, intent(in) :: dims(1:3)
 integer :: coordn(1:3)
 integer :: coordrecv(1:3)
 integer :: ibuf,kbuf,tag,rrank,rankrecv,ierr
 integer :: status(MPI_STATUS_SIZE,2*nproc), status1(MPI_STATUS_SIZE)
 integer :: myrank

 real(8), allocatable, save :: sarr(:,:) ! Send array
 real(8), allocatable, save :: rarr(:,:) ! Receive array
 real(8) :: pr_time, pro_time

 integer, allocatable, save :: maxbuf(:) ! The maximal size of send buffer for each process 
 integer, allocatable, save :: req2(:)

 logical :: firstcall = .true.

! pr_time = MPI_WTIME()
 call cpu_time(pr_time)

 call MPI_COMM_RANK(comm3d, myrank, ierr)

 if (firstcall) then
   allocate (maxbuf(0:nproc-1))
   allocate (req2(1:2*nproc))
   call MAXBUFCALC()
   max_maxbuf = 1
   do rrank = 0, nproc - 1
     max_maxbuf = max(max_maxbuf,maxbuf(rrank))
   enddo
   allocate (sarr(1:maxbuf(myrank), 0:nproc-1) )
   allocate (rarr(1:max_maxbuf,     0:nproc-1) )
 endif

 tag = 1
 coordn(2) = coords(1) ; coordrecv(1) = coords(2)
 coordn(3) = coords(3) ; coordrecv(3) = coords(3)

 ki = nxpi(myrank)
 if (coords(2) /= 0) then
   do i = 0, coords(2)-1
     ki = ki + isize(i)
   enddo
 endif

 do icoord = coords(2), dims(2) - 1
   coordn(1) = icoord
   coordrecv(2) = icoord
   call MPI_CART_RANK(COMM3D,coordn,rrank,ierr)
   call MPI_CART_RANK(COMM3D,coordrecv,rankrecv,ierr)
   ibuf = 1
   do k = nzpi(myrank), nzpe(myrank)
   do j = nypi(myrank), nype(myrank)
   do i = ki, ki + isize(icoord) - 1 ! Sending for coordn(2) = icoord process the isize(icoord) portion
     sarr(ibuf,rrank) = F(i,j,k)
     ibuf = ibuf + 1
   enddo
   enddo
   enddo


   call MPI_ISEND(sarr(1,rrank),maxbuf(myrank),MPI_DOUBLE_PRECISION, &
     &    rrank,tag,comm3d,req2(icoord+nproc+1),ierr)

   call MPI_IRECV(rarr(1,rankrecv),maxbuf(rankrecv),MPI_DOUBLE_PRECISION, &
     &    rankrecv,tag,comm3d,req2(icoord+1),ierr)


   ki = ki + isize(icoord)
 enddo

 ki = nxpi(myrank)
 if (coords(2) /= 0) then
   do icoord = 0, coords(2) - 1
     coordn(1) = icoord
     coordrecv(2) = icoord
     call MPI_CART_RANK(COMM3D,coordn,rrank,ierr)
     call MPI_CART_RANK(COMM3D,coordrecv,rankrecv,ierr)
     ibuf = 1
     do k = nzpi(myrank), nzpe(myrank)
     do j = nypi(myrank), nype(myrank)
     do i = ki, ki + isize(icoord) - 1 ! Sending for coordn(2) = icoord process the isize(icoord) portion
       sarr(ibuf,rrank) = F(i,j,k)
       ibuf = ibuf + 1
     enddo
     enddo
     enddo


     call MPI_ISEND(sarr(1,rrank),maxbuf(myrank),MPI_DOUBLE_PRECISION, &
     &    rrank,tag,comm3d,req2(icoord+nproc+1),ierr)

     call MPI_IRECV(rarr(1,rankrecv),maxbuf(rankrecv),MPI_DOUBLE_PRECISION, &
     &    rankrecv,tag,comm3d,req2(icoord+1),ierr)


    ki = ki + isize(icoord)
  enddo
endif

 call MPI_WAITALL(dims(1), req2(1), status(1,1), ierr)

 coordn(1) = coords(2)
 do icoord = 0, dims(2) - 1
   coordn(2) = icoord
   call MPI_CART_RANK(COMM3D,coordn,rrank,ierr)
   ibuf = 1
   do k = nzpi(rrank), nzpe(rrank)
     do j = nypi(rrank), nype(rrank)
       do i = imin, imin + isize(coords(1)) - 1
         Q(i,j,k) = rarr(ibuf,rrank)
         ibuf = ibuf + 1
       enddo
     enddo
   enddo
 enddo

 call MPI_WAITALL(dims(1), req2(nproc+1), status(1,1), ierr)

 call MPI_BARRIER(comm3d, ierr)

! proc_time(transpose_yz_xz_indic,2) = proc_time(transpose_yz_xz_indic,2) + &
! & MPI_WTIME() - pr_time

 call cpu_time(pro_time)
 proc_time(transpose_yz_xz_indic,2) = proc_time(transpose_yz_xz_indic,2) + &
 & pro_time - pr_time

 if (firstcall) firstcall = .false.
 RETURN
 contains


 SUBROUTINE MAXBUFCALC

 implicit none

 k = 1
 do icoord = 0, dims(1) - 1
   k = max(k, isize(icoord))
 enddo

 do j = 0, dims(2) - 1
   do i = 0, dims(1) - 1
     coordn(1) = i
     coordn(2) = j
     coordn(3) = 0
     call MPI_CART_RANK(COMM3D,coordn,rrank,ierr)
     maxbuf(rrank) = (nzpe(rrank) - nzpi(rrank) + 1) * &
     &               (nype(rrank) - nypi(rrank) + 1) * k
   enddo
 enddo

 END SUBROUTINE MAXBUFCALC
 END SUBROUTINE TRANSPOSE_YZ_XZ

